Inference methods comparison for the simplified Naomi model

Abstract

Background We have run the simplified Naomi model using a range of inference methods: TMB, aghq, adam and tmbstan.

Task In this report, we compare the accuracy of the posterior distributions obtained from these inference methods using histograms, Kolmogorov-Smirnov tests, maximum mean discrepancy and Pareto-smoothed importance sampling.

1 Background

We compare the inference results from TMB, aghq, adam, and tmbstan. Import these inference results as follows:

tmb <- readRDS("depends/tmb.rds")
aghq <- readRDS("depends/aghq.rds")
adam <- readRDS("depends/adam.rds")
tmbstan <- readRDS("depends/tmbstan.rds")

depends <- yaml::read_yaml("orderly.yml")$depends

1.1 Run details

For more information about the conditions under which these results were generated, see:

1.1.1 TMB

dependency_details <- function(i) {
  report_name <- names(depends[[i]])
  print(paste0("Inference results obtained from ", report_name, " with the query ", depends[[i]][[report_name]]$id))
  report_id <- orderly::orderly_search(query = depends[[i]][[report_name]]$id, report_name)
  print(paste0("Obtained report had ID ", report_id, " and was run with the following parameters:"))
  print(orderly::orderly_info(report_id, report_name)$parameters)
}

dependency_details(1)

## [1] "Inference results obtained from naomi-simple_fit with the query latest(parameter:tmb == TRUE)"
## [1] "Obtained report had ID 20230119-142551-23704e7b and was run with the following parameters:"
## $tmb
## [1] TRUE
## 
## $aghq
## [1] FALSE
## 
## $k
## [1] 2
## 
## $ndConstruction
## [1] "product"
## 
## $tmbstan
## [1] FALSE
## 
## $niter
## [1] 1000
## 
## $nthin
## [1] 1
## 
## $nsample
## [1] 1000

1.1.2 aghq

dependency_details(2)

## [1] "Inference results obtained from naomi-simple_fit with the query latest(parameter:aghq == TRUE && parameter:k == 1)"
## [1] "Obtained report had ID 20230119-170459-9b07786e and was run with the following parameters:"
## $aghq
## [1] TRUE
## 
## $k
## [1] 1
## 
## $ndConstruction
## [1] "product"
## 
## $tmb
## [1] FALSE
## 
## $tmbstan
## [1] FALSE
## 
## $niter
## [1] 1000
## 
## $nthin
## [1] 1
## 
## $nsample
## [1] 1000

1.1.3 adam

dependency_details(3)

## [1] "Inference results obtained from naomi-simple_fit with the query latest(parameter:adam == TRUE)"
## [1] "Obtained report had ID 20230209-142804-11535a68 and was run with the following parameters:"
## $adam
## [1] TRUE
## 
## $tmb
## [1] FALSE
## 
## $aghq
## [1] FALSE
## 
## $k
## [1] 2
## 
## $ndConstruction
## [1] "product"
## 
## $tmbstan
## [1] FALSE
## 
## $niter
## [1] 1000
## 
## $nthin
## [1] 1
## 
## $nsample
## [1] 1000
## 
## $area_level
## [1] 4

1.1.4 tmbstan

dependency_details(4)

## [1] "Inference results obtained from naomi-simple_fit with the query latest(parameter:tmbstan == TRUE)"
## [1] "Obtained report had ID 20230114-142916-efabd65d and was run with the following parameters:"
## $tmbstan
## [1] TRUE
## 
## $niter
## [1] 20000
## 
## $nthin
## [1] 20
## 
## $tmb
## [1] FALSE
## 
## $aghq
## [1] FALSE
## 
## $k
## [1] 2
## 
## $nsample
## [1] 1000

1.2 Time taken

data.frame(
  "TMB" = tmb$time,
  "aghq" = aghq$time,
  "adam" = adam$time,
  "tmbstan" = tmbstan$time
)

adam <- adam$adam

2 Histograms and ECDF difference plots

We create histograms and empirical cumulative distribution function (ECDF) difference plots of the samples from each method. All of the possible parameter names are as follows:

pars <- unique(names(tmb$fit$obj$env$par))
pars

##  [1] "beta_rho"             "beta_alpha"           "beta_lambda"          "beta_anc_rho"         "beta_anc_alpha"       "logit_phi_rho_x"      "log_sigma_rho_x"      "logit_phi_rho_xs"    
##  [9] "log_sigma_rho_xs"     "logit_phi_rho_a"      "log_sigma_rho_a"      "logit_phi_rho_as"     "log_sigma_rho_as"     "log_sigma_rho_xa"     "u_rho_x"              "us_rho_x"            
## [17] "u_rho_xs"             "us_rho_xs"            "u_rho_a"              "u_rho_as"             "logit_phi_alpha_x"    "log_sigma_alpha_x"    "logit_phi_alpha_xs"   "log_sigma_alpha_xs"  
## [25] "logit_phi_alpha_a"    "log_sigma_alpha_a"    "logit_phi_alpha_as"   "log_sigma_alpha_as"   "log_sigma_alpha_xa"   "u_alpha_x"            "us_alpha_x"           "u_alpha_xs"          
## [33] "us_alpha_xs"          "u_alpha_a"            "u_alpha_as"           "u_alpha_xa"           "OmegaT_raw"           "log_betaT"            "log_sigma_lambda_x"   "ui_lambda_x"         
## [41] "log_sigma_ancrho_x"   "log_sigma_ancalpha_x" "ui_anc_rho_x"         "ui_anc_alpha_x"       "log_sigma_or_gamma"   "log_or_gamma"

We will produce plots about the following subset of them. There is no particular reason to choose this subset rather than other, it’s quite arbitrary.

pars_eval <- pars %in% c("beta_rho", "beta_alpha", "beta_lambda", "beta_anc_rho", "beta_anc_alpha", "u_rho_x", "us_rho_x", "u_rho_xs")
names(pars_eval) <- pars

beta_rho <- histogram_and_ecdf("beta_rho", i = 1, return_df = TRUE)
saveRDS(beta_rho$df, "beta_rho.rds")

2.1 beta_rho

histogram_and_ecdf_helper <- function(par) lapply(1:sum(names(tmb$fit$obj$env$par) == par), histogram_and_ecdf, par = par)
histogram_and_ecdf_helper("beta_rho")

## [[1]]

## 
## [[2]]

2.2 beta_alpha

histogram_and_ecdf_helper("beta_alpha")

## [[1]]

## 
## [[2]]

2.3 beta_lambda

histogram_and_ecdf_helper("beta_lambda")

## [[1]]

## 
## [[2]]

2.4 beta_anc_rho

histogram_and_ecdf("beta_anc_rho")

2.5 beta_anc_alpha

histogram_and_ecdf("beta_anc_alpha")

2.6 logit

lapply(pars[stringr::str_starts(pars, "logit")], histogram_and_ecdf)

2.7 log_sigma

lapply(pars[stringr::str_starts(pars, "log_sigma")], histogram_and_ecdf)

2.8 u_rho_x

histogram_and_ecdf_helper("u_rho_x")

## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

## 
## [[19]]

## 
## [[20]]

## 
## [[21]]

## 
## [[22]]

## 
## [[23]]

## 
## [[24]]

## 
## [[25]]

## 
## [[26]]

## 
## [[27]]

## 
## [[28]]

## 
## [[29]]

## 
## [[30]]

## 
## [[31]]

## 
## [[32]]

2.9 u_rho_x

histogram_and_ecdf_helper("u_rho_x")

## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

## 
## [[19]]

## 
## [[20]]

## 
## [[21]]

## 
## [[22]]

## 
## [[23]]

## 
## [[24]]

## 
## [[25]]

## 
## [[26]]

## 
## [[27]]

## 
## [[28]]

## 
## [[29]]

## 
## [[30]]

## 
## [[31]]

## 
## [[32]]

2.10 us_rho_x

histogram_and_ecdf_helper("us_rho_x")

## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

## 
## [[19]]

## 
## [[20]]

## 
## [[21]]

## 
## [[22]]

## 
## [[23]]

## 
## [[24]]

## 
## [[25]]

## 
## [[26]]

## 
## [[27]]

## 
## [[28]]

## 
## [[29]]

## 
## [[30]]

## 
## [[31]]

## 
## [[32]]

2.11 u_rho_xs

histogram_and_ecdf_helper("u_rho_xs")

## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

## 
## [[19]]

## 
## [[20]]

## 
## [[21]]

## 
## [[22]]

## 
## [[23]]

## 
## [[24]]

## 
## [[25]]

## 
## [[26]]

## 
## [[27]]

## 
## [[28]]

## 
## [[29]]

## 
## [[30]]

## 
## [[31]]

## 
## [[32]]

2.12 us_rho_xs

histogram_and_ecdf_helper("us_rho_xs")

2.13 u_rho_as

histogram_and_ecdf_helper("u_rho_as")

2.14 u_alpha_x

histogram_and_ecdf_helper("u_alpha_x")

2.15 us_alpha_x

histogram_and_ecdf_helper("us_alpha_x")

2.16 u_alpha_xs

histogram_and_ecdf_helper("u_alpha_xs")

2.17 us_alpha_xs

histogram_and_ecdf_helper("us_alpha_xs")

## [[1]]

## 
## [[2]]

## 
## [[3]]

## 
## [[4]]

## 
## [[5]]

## 
## [[6]]

## 
## [[7]]

## 
## [[8]]

## 
## [[9]]

## 
## [[10]]

## 
## [[11]]

## 
## [[12]]

## 
## [[13]]

## 
## [[14]]

## 
## [[15]]

## 
## [[16]]

## 
## [[17]]

## 
## [[18]]

## 
## [[19]]

## 
## [[20]]

## 
## [[21]]

## 
## [[22]]

## 
## [[23]]

## 
## [[24]]

## 
## [[25]]

## 
## [[26]]

## 
## [[27]]

## 
## [[28]]

## 
## [[29]]

## 
## [[30]]

## 
## [[31]]

## 
## [[32]]

2.18 u_alpha_a

histogram_and_ecdf_helper("u_alpha_a")

2.19 u_alpha_as

histogram_and_ecdf_helper("u_alpha_as")

2.20 u_alpha_xa

histogram_and_ecdf_helper("u_alpha_xa")

2.21 ui_lambda_x

histogram_and_ecdf_helper("ui_lambda_x")

2.22 ui_anc_rho_x

histogram_and_ecdf_helper("ui_anc_rho_x")

2.23 ui_anc_alpha_x

histogram_and_ecdf_helper("ui_anc_alpha_x")

2.24 log_or_gamma

histogram_and_ecdf_helper("log_or_gamma")

3 KS plots

3.1 Individual parameters

ks_helper <- function(par, sw = FALSE, ...) {
  to_ks_df(par, starts_with = sw) %>% ks_plot(par, ...)
}

3.1.1 beta

ks_helper("beta", sw = TRUE)

ks_helper("beta", sw = TRUE, method1 = "TMB", method2 = "TMB")

3.1.2 logit

ks_helper("logit", sw = TRUE)

ks_helper("logit", sw = TRUE, method1 = "TMB", method2 = "adam")

3.1.3 log_sigma

ks_helper("log_sigma", sw = TRUE)

ks_helper("log_sigma", sw = TRUE , method1 = "TMB", method2 = "adam")

3.1.4 u_rho_x

ks_helper("u_rho_x")

## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

ks_helper("u_rho_x", method1 = "TMB", method2 = "adam")

## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

3.1.5 u_rho_xs

ks_helper("u_rho_xs")

## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

ks_helper("u_rho_xs", method1 = "TMB", method2 = "adam")

## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

3.1.6 us_rho_x

ks_helper("us_rho_x")

## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

ks_helper("us_rho_x", method1 = "TMB", method2 = "adam")

## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`

3.1.7 us_rho_xs

ks_helper("us_rho_xs")
ks_helper("us_rho_xs", method1 = "TMB", method2 = "adam")

3.1.8 u_rho_a

ks_helper("u_rho_a")
ks_helper("u_rho_a", method1 = "TMB", method2 = "adam")

3.1.9 u_rho_as

ks_helper("u_rho_as")
ks_helper("u_rho_as", method1 = "TMB", method2 = "adam")

3.1.10 u_alpha_x

ks_helper("u_alpha_x")
ks_helper("u_alpha_x", method1 = "TMB", method2 = "adam")

3.1.11 u_alpha_xs

ks_helper("u_alpha_xs")
ks_helper("u_alpha_xs", method1 = "TMB", method2 = "adam")

3.1.12 us_alpha_x

ks_helper("us_alpha_x")
ks_helper("us_alpha_x", method1 = "TMB", method2 = "adam")

3.1.13 us_alpha_xs

ks_helper("us_alpha_xs")
ks_helper("us_alpha_xs", method1 = "TMB", method2 = "adam")

3.1.14 u_alpha_a

ks_helper("u_alpha_a")
ks_helper("u_alpha_a", method1 = "TMB", method2 = "adam")

3.1.15 u_alpha_as

ks_helper("u_alpha_as")
ks_helper("u_alpha_as", method1 = "TMB", method2 = "adam")

3.1.16 u_alpha_xa

ks_helper("u_alpha_xa")
ks_helper("u_alpha_xa", method1 = "TMB", method2 = "adam")

3.1.17 ui_anc_rho_x

ks_helper("ui_anc_rho_x")
ks_helper("ui_anc_rho_x", method1 = "TMB", method2 = "adam")

3.1.18 ui_anc_alpha_x

ks_helper("ui_anc_alpha_x")
ks_helper("ui_anc_alpha_x", method1 = "TMB", method2 = "adam")

3.1.19 log_or_gamma

ks_helper("log_or_gamma")
ks_helper("log_or_gamma", method1 = "TMB", method2 = "adam")

3.2 Summary table

options(dplyr.summarise.inform = FALSE)

ks_summary <- lapply(unique(names(tmb$fit$obj$env$par)), function(x) {
  to_ks_df(x) %>%
    group_by(method, index) %>%
    summarise(ks = mean(ks), par = x) %>%
    ungroup()
}) %>%
  bind_rows() %>%
  pivot_wider(names_from = "method", values_from = "ks") %>%
  rename(
    "Parameter" = "par",
    "KS(adam, tmbstan)" = "adam",
    "KS(aghq, tmbstan)" = "aghq",
    "KS(TMB, tmbstan)" = "TMB",
  )

## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## New names:
## • `` -> `...2`

r <- adam$quad$obj$env$random
x_names <- names(adam$quad$obj$env$par[r])
theta_names <- names(adam$quad$obj$env$par[-r])

dict <- data.frame(
  Parameter = c(unique(x_names), unique(theta_names)),
  Type = c(rep("Latent field", length(unique(x_names))), rep("Hyper", length(unique(theta_names))))
)

ks_summary <- ks_summary %>%
  left_join(dict, by = "Parameter")

Values in green are the lowest KS difference for that particular parameter:

ks_summary %>%
  gt::gt() %>%
  gt::fmt_number(
    columns = starts_with("KS"),
    decimals = 3
  ) %>%
  gt::tab_style(
    style = cell_fill(color = "#009E73"),
    locations = cells_body(
      columns = `KS(adam, tmbstan)`,
      rows = `KS(adam, tmbstan)` < `KS(aghq, tmbstan)` & `KS(adam, tmbstan)` < `KS(TMB, tmbstan)`
    )
  ) %>%
    gt::tab_style(
    style = cell_fill(color = "#009E73"),
    locations = cells_body(
      columns = `KS(aghq, tmbstan)`,
      rows = `KS(aghq, tmbstan)` < `KS(adam, tmbstan)` & `KS(aghq, tmbstan)` < `KS(TMB, tmbstan)`
    )
  ) %>%
    gt::tab_style(
    style = cell_fill(color = "#009E73"),
    locations = cells_body(
      columns = `KS(TMB, tmbstan)`,
      rows = `KS(TMB, tmbstan)` < `KS(aghq, tmbstan)` & `KS(TMB, tmbstan)` < `KS(adam, tmbstan)`
    )
  )

index	Parameter	KS(adam, tmbstan)	KS(aghq, tmbstan)	KS(TMB, tmbstan)	Type
1	beta_rho	0.104	0.091	0.092	Latent field
2	beta_rho	0.079	0.088	0.078	Latent field
1	beta_alpha	0.094	0.084	0.082	Latent field
2	beta_alpha	0.152	0.118	0.124	Latent field
1	beta_lambda	0.062	0.096	0.113	Latent field
2	beta_lambda	0.022	0.027	0.036	Latent field
1	beta_anc_rho	0.069	0.105	0.090	Latent field
1	beta_anc_alpha	0.070	0.044	0.025	Latent field
1	logit_phi_rho_x	0.562	0.562	0.125	Hyper
1	log_sigma_rho_x	0.656	0.656	0.158	Hyper
1	logit_phi_rho_xs	0.510	0.510	0.111	Hyper
1	log_sigma_rho_xs	0.613	0.613	0.203	Hyper
1	logit_phi_rho_a	0.536	0.536	0.056	Hyper
1	log_sigma_rho_a	0.581	0.581	0.092	Hyper
1	logit_phi_rho_as	0.569	0.569	0.104	Hyper
1	log_sigma_rho_as	0.589	0.589	0.134	Hyper
1	log_sigma_rho_xa	0.669	0.669	0.224	Hyper
1	u_rho_x	0.040	0.096	0.065	Latent field
2	u_rho_x	0.019	0.037	0.052	Latent field
3	u_rho_x	0.042	0.060	0.067	Latent field
4	u_rho_x	0.046	0.077	0.082	Latent field
5	u_rho_x	0.075	0.079	0.069	Latent field
6	u_rho_x	0.112	0.091	0.104	Latent field
7	u_rho_x	0.044	0.088	0.108	Latent field
8	u_rho_x	0.070	0.052	0.064	Latent field
9	u_rho_x	0.058	0.096	0.089	Latent field
10	u_rho_x	0.040	0.080	0.054	Latent field
11	u_rho_x	0.053	0.044	0.064	Latent field
12	u_rho_x	0.059	0.132	0.116	Latent field
13	u_rho_x	0.034	0.056	0.031	Latent field
14	u_rho_x	0.081	0.092	0.103	Latent field
15	u_rho_x	0.084	0.040	0.055	Latent field
16	u_rho_x	0.026	0.022	0.021	Latent field
17	u_rho_x	0.035	0.041	0.062	Latent field
18	u_rho_x	0.059	0.083	0.061	Latent field
19	u_rho_x	0.024	0.026	0.071	Latent field
20	u_rho_x	0.022	0.045	0.040	Latent field
21	u_rho_x	0.036	0.052	0.044	Latent field
22	u_rho_x	0.119	0.107	0.089	Latent field
23	u_rho_x	0.043	0.055	0.043	Latent field
24	u_rho_x	0.064	0.033	0.059	Latent field
25	u_rho_x	0.063	0.087	0.078	Latent field
26	u_rho_x	0.065	0.025	0.030	Latent field
27	u_rho_x	0.060	0.083	0.100	Latent field
28	u_rho_x	0.053	0.050	0.073	Latent field
29	u_rho_x	0.088	0.106	0.126	Latent field
30	u_rho_x	0.078	0.060	0.066	Latent field
31	u_rho_x	0.051	0.042	0.038	Latent field
32	u_rho_x	0.037	0.088	0.068	Latent field
1	us_rho_x	0.162	0.165	0.139	Latent field
2	us_rho_x	0.166	0.139	0.123	Latent field
3	us_rho_x	0.088	0.110	0.091	Latent field
4	us_rho_x	0.075	0.092	0.092	Latent field
5	us_rho_x	0.073	0.087	0.082	Latent field
6	us_rho_x	0.188	0.190	0.103	Latent field
7	us_rho_x	0.094	0.057	0.081	Latent field
8	us_rho_x	0.083	0.070	0.073	Latent field
9	us_rho_x	0.046	0.091	0.057	Latent field
10	us_rho_x	0.066	0.042	0.059	Latent field
11	us_rho_x	0.104	0.101	0.091	Latent field
12	us_rho_x	0.073	0.086	0.068	Latent field
13	us_rho_x	0.076	0.073	0.065	Latent field
14	us_rho_x	0.066	0.055	0.029	Latent field
15	us_rho_x	0.043	0.077	0.052	Latent field
16	us_rho_x	0.039	0.075	0.053	Latent field
17	us_rho_x	0.068	0.083	0.082	Latent field
18	us_rho_x	0.115	0.126	0.106	Latent field
19	us_rho_x	0.073	0.095	0.101	Latent field
20	us_rho_x	0.133	0.166	0.124	Latent field
21	us_rho_x	0.102	0.105	0.120	Latent field
22	us_rho_x	0.072	0.063	0.051	Latent field
23	us_rho_x	0.107	0.071	0.080	Latent field
24	us_rho_x	0.074	0.060	0.094	Latent field
25	us_rho_x	0.071	0.084	0.112	Latent field
26	us_rho_x	0.101	0.091	0.074	Latent field
27	us_rho_x	0.042	0.037	0.051	Latent field
28	us_rho_x	0.074	0.043	0.077	Latent field
29	us_rho_x	0.071	0.076	0.088	Latent field
30	us_rho_x	0.084	0.083	0.070	Latent field
31	us_rho_x	0.038	0.033	0.044	Latent field
32	us_rho_x	0.090	0.070	0.073	Latent field
1	u_rho_xs	0.081	0.082	0.093	Latent field
2	u_rho_xs	0.095	0.081	0.104	Latent field
3	u_rho_xs	0.081	0.089	0.088	Latent field
4	u_rho_xs	0.088	0.069	0.086	Latent field
5	u_rho_xs	0.169	0.153	0.142	Latent field
6	u_rho_xs	0.187	0.199	0.182	Latent field
7	u_rho_xs	0.104	0.101	0.113	Latent field
8	u_rho_xs	0.110	0.079	0.094	Latent field
9	u_rho_xs	0.118	0.143	0.147	Latent field
10	u_rho_xs	0.088	0.108	0.141	Latent field
11	u_rho_xs	0.131	0.087	0.084	Latent field
12	u_rho_xs	0.189	0.186	0.152	Latent field
13	u_rho_xs	0.120	0.090	0.109	Latent field
14	u_rho_xs	0.136	0.138	0.131	Latent field
15	u_rho_xs	0.127	0.154	0.166	Latent field
16	u_rho_xs	0.099	0.076	0.092	Latent field
17	u_rho_xs	0.109	0.081	0.095	Latent field
18	u_rho_xs	0.136	0.125	0.120	Latent field
19	u_rho_xs	0.121	0.133	0.113	Latent field
20	u_rho_xs	0.096	0.104	0.102	Latent field
21	u_rho_xs	0.141	0.091	0.103	Latent field
22	u_rho_xs	0.213	0.194	0.154	Latent field
23	u_rho_xs	0.101	0.114	0.104	Latent field
24	u_rho_xs	0.137	0.122	0.139	Latent field
25	u_rho_xs	0.093	0.079	0.084	Latent field
26	u_rho_xs	0.117	0.104	0.153	Latent field
27	u_rho_xs	0.149	0.165	0.145	Latent field
28	u_rho_xs	0.142	0.133	0.120	Latent field
29	u_rho_xs	0.201	0.192	0.191	Latent field
30	u_rho_xs	0.135	0.111	0.116	Latent field
31	u_rho_xs	0.067	0.082	0.112	Latent field
32	u_rho_xs	0.084	0.093	0.104	Latent field
1	us_rho_xs	0.046	0.025	0.040	Latent field
2	us_rho_xs	0.048	0.023	0.031	Latent field
3	us_rho_xs	0.037	0.034	0.035	Latent field
4	us_rho_xs	0.031	0.035	0.035	Latent field
5	us_rho_xs	0.036	0.043	0.036	Latent field
6	us_rho_xs	0.025	0.031	0.022	Latent field
7	us_rho_xs	0.025	0.046	0.018	Latent field
8	us_rho_xs	0.048	0.044	0.043	Latent field
9	us_rho_xs	0.054	0.037	0.060	Latent field
10	us_rho_xs	0.042	0.052	0.058	Latent field
11	us_rho_xs	0.041	0.042	0.061	Latent field
12	us_rho_xs	0.043	0.036	0.065	Latent field
13	us_rho_xs	0.043	0.042	0.061	Latent field
14	us_rho_xs	0.052	0.033	0.061	Latent field
15	us_rho_xs	0.040	0.024	0.053	Latent field
16	us_rho_xs	0.062	0.022	0.041	Latent field
17	us_rho_xs	0.056	0.045	0.038	Latent field
18	us_rho_xs	0.051	0.036	0.040	Latent field
19	us_rho_xs	0.035	0.029	0.040	Latent field
20	us_rho_xs	0.059	0.042	0.057	Latent field
21	us_rho_xs	0.065	0.038	0.049	Latent field
22	us_rho_xs	0.047	0.051	0.048	Latent field
23	us_rho_xs	0.034	0.040	0.043	Latent field
24	us_rho_xs	0.034	0.028	0.026	Latent field
25	us_rho_xs	0.032	0.042	0.044	Latent field
26	us_rho_xs	0.044	0.029	0.047	Latent field
27	us_rho_xs	0.051	0.049	0.040	Latent field
28	us_rho_xs	0.022	0.041	0.034	Latent field
29	us_rho_xs	0.062	0.031	0.038	Latent field
30	us_rho_xs	0.055	0.033	0.032	Latent field
31	us_rho_xs	0.073	0.036	0.049	Latent field
32	us_rho_xs	0.036	0.033	0.052	Latent field
1	u_rho_a	0.085	0.091	0.092	Latent field
2	u_rho_a	0.078	0.089	0.088	Latent field
3	u_rho_a	0.089	0.091	0.092	Latent field
4	u_rho_a	0.095	0.090	0.089	Latent field
5	u_rho_a	0.088	0.090	0.090	Latent field
6	u_rho_a	0.091	0.091	0.091	Latent field
7	u_rho_a	0.082	0.089	0.091	Latent field
8	u_rho_a	0.088	0.087	0.089	Latent field
9	u_rho_a	0.097	0.088	0.090	Latent field
10	u_rho_a	0.089	0.084	0.087	Latent field
1	u_rho_as	0.096	0.082	0.071	Latent field
2	u_rho_as	0.081	0.082	0.072	Latent field
3	u_rho_as	0.102	0.084	0.076	Latent field
4	u_rho_as	0.092	0.086	0.075	Latent field
5	u_rho_as	0.091	0.086	0.079	Latent field
6	u_rho_as	0.094	0.083	0.077	Latent field
7	u_rho_as	0.096	0.087	0.080	Latent field
8	u_rho_as	0.085	0.088	0.079	Latent field
9	u_rho_as	0.091	0.086	0.076	Latent field
10	u_rho_as	0.090	0.082	0.076	Latent field
1	logit_phi_alpha_x	0.512	0.512	0.071	Hyper
1	log_sigma_alpha_x	0.505	0.505	0.125	Hyper
1	logit_phi_alpha_xs	0.551	0.551	0.143	Hyper
1	log_sigma_alpha_xs	0.540	0.540	0.159	Hyper
1	logit_phi_alpha_a	0.533	0.533	0.085	Hyper
1	log_sigma_alpha_a	0.569	0.569	0.124	Hyper
1	logit_phi_alpha_as	0.516	0.516	0.081	Hyper
1	log_sigma_alpha_as	0.514	0.514	0.098	Hyper
1	log_sigma_alpha_xa	0.627	0.627	0.146	Hyper
1	u_alpha_x	0.089	0.149	0.092	Latent field
2	u_alpha_x	0.135	0.134	0.079	Latent field
3	u_alpha_x	0.120	0.115	0.128	Latent field
4	u_alpha_x	0.104	0.144	0.129	Latent field
5	u_alpha_x	0.151	0.125	0.129	Latent field
6	u_alpha_x	0.095	0.085	0.101	Latent field
7	u_alpha_x	0.042	0.067	0.064	Latent field
8	u_alpha_x	0.130	0.118	0.115	Latent field
9	u_alpha_x	0.146	0.156	0.166	Latent field
10	u_alpha_x	0.171	0.142	0.095	Latent field
11	u_alpha_x	0.132	0.144	0.153	Latent field
12	u_alpha_x	0.127	0.076	0.088	Latent field
13	u_alpha_x	0.124	0.138	0.125	Latent field
14	u_alpha_x	0.081	0.090	0.112	Latent field
15	u_alpha_x	0.090	0.096	0.085	Latent field
16	u_alpha_x	0.077	0.064	0.087	Latent field
17	u_alpha_x	0.054	0.074	0.058	Latent field
18	u_alpha_x	0.102	0.115	0.082	Latent field
19	u_alpha_x	0.104	0.139	0.130	Latent field
20	u_alpha_x	0.101	0.111	0.106	Latent field
21	u_alpha_x	0.170	0.155	0.100	Latent field
22	u_alpha_x	0.125	0.141	0.156	Latent field
23	u_alpha_x	0.153	0.176	0.122	Latent field
24	u_alpha_x	0.138	0.156	0.168	Latent field
25	u_alpha_x	0.082	0.114	0.117	Latent field
26	u_alpha_x	0.096	0.066	0.069	Latent field
27	u_alpha_x	0.115	0.150	0.142	Latent field
28	u_alpha_x	0.121	0.108	0.101	Latent field
29	u_alpha_x	0.102	0.113	0.150	Latent field
30	u_alpha_x	0.103	0.132	0.132	Latent field
31	u_alpha_x	0.093	0.104	0.089	Latent field
32	u_alpha_x	0.064	0.069	0.097	Latent field
1	us_alpha_x	0.063	0.094	0.053	Latent field
2	us_alpha_x	0.073	0.077	0.048	Latent field
3	us_alpha_x	0.057	0.064	0.096	Latent field
4	us_alpha_x	0.050	0.044	0.066	Latent field
5	us_alpha_x	0.061	0.083	0.075	Latent field
6	us_alpha_x	0.096	0.083	0.064	Latent field
7	us_alpha_x	0.036	0.023	0.038	Latent field
8	us_alpha_x	0.092	0.087	0.087	Latent field
9	us_alpha_x	0.117	0.090	0.118	Latent field
10	us_alpha_x	0.107	0.107	0.102	Latent field
11	us_alpha_x	0.118	0.067	0.092	Latent field
12	us_alpha_x	0.056	0.028	0.063	Latent field
13	us_alpha_x	0.088	0.102	0.098	Latent field
14	us_alpha_x	0.032	0.029	0.040	Latent field
15	us_alpha_x	0.104	0.049	0.067	Latent field
16	us_alpha_x	0.030	0.048	0.024	Latent field
17	us_alpha_x	0.057	0.083	0.081	Latent field
18	us_alpha_x	0.044	0.072	0.089	Latent field
19	us_alpha_x	0.069	0.087	0.057	Latent field
20	us_alpha_x	0.056	0.059	0.074	Latent field
21	us_alpha_x	0.106	0.106	0.094	Latent field
22	us_alpha_x	0.086	0.127	0.131	Latent field
23	us_alpha_x	0.118	0.101	0.083	Latent field
24	us_alpha_x	0.104	0.087	0.108	Latent field
25	us_alpha_x	0.055	0.083	0.098	Latent field
26	us_alpha_x	0.092	0.089	0.064	Latent field
27	us_alpha_x	0.062	0.054	0.075	Latent field
28	us_alpha_x	0.042	0.060	0.035	Latent field
29	us_alpha_x	0.066	0.061	0.072	Latent field
30	us_alpha_x	0.052	0.069	0.056	Latent field
31	us_alpha_x	0.037	0.040	0.042	Latent field
32	us_alpha_x	0.065	0.060	0.029	Latent field
1	u_alpha_xs	0.074	0.074	0.063	Latent field
2	u_alpha_xs	0.073	0.105	0.060	Latent field
3	u_alpha_xs	0.059	0.075	0.068	Latent field
4	u_alpha_xs	0.067	0.091	0.133	Latent field
5	u_alpha_xs	0.132	0.129	0.124	Latent field
6	u_alpha_xs	0.059	0.065	0.090	Latent field
7	u_alpha_xs	0.083	0.062	0.058	Latent field
8	u_alpha_xs	0.093	0.058	0.050	Latent field
9	u_alpha_xs	0.089	0.090	0.087	Latent field
10	u_alpha_xs	0.123	0.138	0.154	Latent field
11	u_alpha_xs	0.140	0.109	0.092	Latent field
12	u_alpha_xs	0.093	0.117	0.106	Latent field
13	u_alpha_xs	0.105	0.104	0.086	Latent field
14	u_alpha_xs	0.060	0.041	0.052	Latent field
15	u_alpha_xs	0.135	0.095	0.135	Latent field
16	u_alpha_xs	0.167	0.152	0.081	Latent field
17	u_alpha_xs	0.058	0.037	0.067	Latent field
18	u_alpha_xs	0.178	0.168	0.088	Latent field
19	u_alpha_xs	0.081	0.062	0.065	Latent field
20	u_alpha_xs	0.061	0.064	0.083	Latent field
21	u_alpha_xs	0.128	0.105	0.085	Latent field
22	u_alpha_xs	0.193	0.180	0.134	Latent field
23	u_alpha_xs	0.100	0.093	0.116	Latent field
24	u_alpha_xs	0.152	0.132	0.145	Latent field
25	u_alpha_xs	0.184	0.134	0.132	Latent field
26	u_alpha_xs	0.089	0.052	0.083	Latent field
27	u_alpha_xs	0.133	0.128	0.148	Latent field
28	u_alpha_xs	0.107	0.106	0.101	Latent field
29	u_alpha_xs	0.190	0.156	0.106	Latent field
30	u_alpha_xs	0.090	0.089	0.083	Latent field
31	u_alpha_xs	0.105	0.097	0.110	Latent field
32	u_alpha_xs	0.125	0.110	0.090	Latent field
1	us_alpha_xs	0.053	0.044	0.030	Latent field
2	us_alpha_xs	0.025	0.032	0.022	Latent field
3	us_alpha_xs	0.041	0.042	0.072	Latent field
4	us_alpha_xs	0.049	0.064	0.100	Latent field
5	us_alpha_xs	0.089	0.085	0.120	Latent field
6	us_alpha_xs	0.055	0.075	0.086	Latent field
7	us_alpha_xs	0.041	0.018	0.036	Latent field
8	us_alpha_xs	0.096	0.124	0.118	Latent field
9	us_alpha_xs	0.129	0.120	0.119	Latent field
10	us_alpha_xs	0.130	0.149	0.108	Latent field
11	us_alpha_xs	0.163	0.121	0.124	Latent field
12	us_alpha_xs	0.116	0.142	0.135	Latent field
13	us_alpha_xs	0.082	0.096	0.097	Latent field
14	us_alpha_xs	0.106	0.097	0.101	Latent field
15	us_alpha_xs	0.133	0.116	0.123	Latent field
16	us_alpha_xs	0.146	0.126	0.091	Latent field
17	us_alpha_xs	0.063	0.060	0.062	Latent field
18	us_alpha_xs	0.104	0.072	0.059	Latent field
19	us_alpha_xs	0.034	0.025	0.022	Latent field
20	us_alpha_xs	0.046	0.040	0.036	Latent field
21	us_alpha_xs	0.067	0.065	0.054	Latent field
22	us_alpha_xs	0.140	0.113	0.097	Latent field
23	us_alpha_xs	0.113	0.103	0.116	Latent field
24	us_alpha_xs	0.173	0.153	0.145	Latent field
25	us_alpha_xs	0.095	0.115	0.088	Latent field
26	us_alpha_xs	0.144	0.131	0.112	Latent field
27	us_alpha_xs	0.149	0.138	0.144	Latent field
28	us_alpha_xs	0.141	0.129	0.092	Latent field
29	us_alpha_xs	0.153	0.153	0.132	Latent field
30	us_alpha_xs	0.167	0.126	0.153	Latent field
31	us_alpha_xs	0.114	0.110	0.129	Latent field
32	us_alpha_xs	0.113	0.110	0.105	Latent field
1	u_alpha_a	0.052	0.065	0.061	Latent field
2	u_alpha_a	0.054	0.066	0.065	Latent field
3	u_alpha_a	0.068	0.078	0.066	Latent field
4	u_alpha_a	0.062	0.079	0.066	Latent field
5	u_alpha_a	0.056	0.082	0.067	Latent field
6	u_alpha_a	0.058	0.071	0.063	Latent field
7	u_alpha_a	0.068	0.067	0.060	Latent field
8	u_alpha_a	0.054	0.071	0.066	Latent field
9	u_alpha_a	0.053	0.072	0.077	Latent field
10	u_alpha_a	0.056	0.087	0.075	Latent field
11	u_alpha_a	0.075	0.092	0.084	Latent field
12	u_alpha_a	0.077	0.095	0.081	Latent field
13	u_alpha_a	0.079	0.067	0.074	Latent field
1	u_alpha_as	0.082	0.107	0.098	Latent field
2	u_alpha_as	0.109	0.127	0.117	Latent field
3	u_alpha_as	0.140	0.128	0.138	Latent field
4	u_alpha_as	0.162	0.132	0.129	Latent field
5	u_alpha_as	0.151	0.122	0.129	Latent field
6	u_alpha_as	0.132	0.098	0.114	Latent field
7	u_alpha_as	0.127	0.085	0.099	Latent field
8	u_alpha_as	0.075	0.060	0.091	Latent field
9	u_alpha_as	0.088	0.051	0.057	Latent field
10	u_alpha_as	0.051	0.051	0.047	Latent field
1	u_alpha_xa	0.028	0.033	0.053	Latent field
2	u_alpha_xa	0.037	0.043	0.025	Latent field
3	u_alpha_xa	0.033	0.028	0.028	Latent field
4	u_alpha_xa	0.058	0.075	0.062	Latent field
5	u_alpha_xa	0.104	0.123	0.095	Latent field
6	u_alpha_xa	0.082	0.064	0.079	Latent field
7	u_alpha_xa	0.033	0.039	0.038	Latent field
8	u_alpha_xa	0.069	0.052	0.052	Latent field
9	u_alpha_xa	0.065	0.044	0.044	Latent field
10	u_alpha_xa	0.082	0.116	0.101	Latent field
11	u_alpha_xa	0.107	0.113	0.089	Latent field
12	u_alpha_xa	0.117	0.114	0.095	Latent field
13	u_alpha_xa	0.062	0.085	0.050	Latent field
14	u_alpha_xa	0.058	0.081	0.045	Latent field
15	u_alpha_xa	0.057	0.036	0.045	Latent field
16	u_alpha_xa	0.096	0.098	0.061	Latent field
17	u_alpha_xa	0.025	0.033	0.039	Latent field
18	u_alpha_xa	0.078	0.054	0.029	Latent field
19	u_alpha_xa	0.071	0.078	0.071	Latent field
20	u_alpha_xa	0.032	0.058	0.041	Latent field
21	u_alpha_xa	0.088	0.071	0.036	Latent field
22	u_alpha_xa	0.201	0.197	0.181	Latent field
23	u_alpha_xa	0.027	0.048	0.050	Latent field
24	u_alpha_xa	0.077	0.068	0.088	Latent field
25	u_alpha_xa	0.104	0.071	0.069	Latent field
26	u_alpha_xa	0.058	0.048	0.051	Latent field
27	u_alpha_xa	0.086	0.123	0.086	Latent field
28	u_alpha_xa	0.053	0.043	0.030	Latent field
29	u_alpha_xa	0.060	0.037	0.040	Latent field
30	u_alpha_xa	0.060	0.074	0.054	Latent field
31	u_alpha_xa	0.034	0.032	0.033	Latent field
32	u_alpha_xa	0.035	0.045	0.026	Latent field
1	OmegaT_raw	0.518	0.518	0.022	Hyper
1	log_betaT	0.689	0.689	0.244	Hyper
1	log_sigma_lambda_x	0.736	0.736	0.269	Hyper
1	ui_lambda_x	0.179	0.151	0.170	Latent field
2	ui_lambda_x	0.147	0.132	0.132	Latent field
3	ui_lambda_x	0.169	0.170	0.139	Latent field
4	ui_lambda_x	0.187	0.168	0.178	Latent field
5	ui_lambda_x	0.165	0.167	0.137	Latent field
6	ui_lambda_x	0.157	0.132	0.144	Latent field
7	ui_lambda_x	0.154	0.145	0.151	Latent field
8	ui_lambda_x	0.138	0.143	0.166	Latent field
9	ui_lambda_x	0.178	0.170	0.165	Latent field
10	ui_lambda_x	0.202	0.247	0.243	Latent field
11	ui_lambda_x	0.146	0.146	0.141	Latent field
12	ui_lambda_x	0.161	0.171	0.173	Latent field
13	ui_lambda_x	0.182	0.208	0.191	Latent field
14	ui_lambda_x	0.135	0.135	0.146	Latent field
15	ui_lambda_x	0.167	0.208	0.203	Latent field
16	ui_lambda_x	0.188	0.143	0.178	Latent field
17	ui_lambda_x	0.192	0.180	0.189	Latent field
18	ui_lambda_x	0.168	0.149	0.159	Latent field
19	ui_lambda_x	0.177	0.175	0.174	Latent field
20	ui_lambda_x	0.138	0.163	0.157	Latent field
21	ui_lambda_x	0.143	0.151	0.166	Latent field
22	ui_lambda_x	0.143	0.133	0.126	Latent field
23	ui_lambda_x	0.194	0.164	0.184	Latent field
24	ui_lambda_x	0.166	0.145	0.182	Latent field
25	ui_lambda_x	0.144	0.173	0.134	Latent field
26	ui_lambda_x	0.184	0.152	0.186	Latent field
27	ui_lambda_x	0.149	0.125	0.137	Latent field
28	ui_lambda_x	0.140	0.134	0.139	Latent field
29	ui_lambda_x	0.152	0.154	0.140	Latent field
30	ui_lambda_x	0.141	0.161	0.121	Latent field
31	ui_lambda_x	0.165	0.133	0.139	Latent field
32	ui_lambda_x	0.167	0.138	0.143	Latent field
1	log_sigma_ancrho_x	0.504	0.504	0.037	Hyper
1	log_sigma_ancalpha_x	0.519	0.519	0.083	Hyper
1	ui_anc_rho_x	0.040	0.085	0.072	Latent field
2	ui_anc_rho_x	0.030	0.033	0.019	Latent field
3	ui_anc_rho_x	0.044	0.038	0.090	Latent field
4	ui_anc_rho_x	0.041	0.104	0.092	Latent field
5	ui_anc_rho_x	0.022	0.064	0.060	Latent field
6	ui_anc_rho_x	0.104	0.053	0.053	Latent field
7	ui_anc_rho_x	0.040	0.082	0.075	Latent field
8	ui_anc_rho_x	0.099	0.040	0.051	Latent field
9	ui_anc_rho_x	0.030	0.042	0.040	Latent field
10	ui_anc_rho_x	0.029	0.125	0.125	Latent field
11	ui_anc_rho_x	0.023	0.025	0.055	Latent field
12	ui_anc_rho_x	0.051	0.064	0.037	Latent field
13	ui_anc_rho_x	0.026	0.036	0.030	Latent field
14	ui_anc_rho_x	0.040	0.050	0.067	Latent field
15	ui_anc_rho_x	0.028	0.076	0.044	Latent field
16	ui_anc_rho_x	0.041	0.027	0.026	Latent field
17	ui_anc_rho_x	0.038	0.042	0.057	Latent field
18	ui_anc_rho_x	0.063	0.032	0.033	Latent field
19	ui_anc_rho_x	0.037	0.042	0.082	Latent field
20	ui_anc_rho_x	0.026	0.038	0.030	Latent field
21	ui_anc_rho_x	0.047	0.046	0.065	Latent field
22	ui_anc_rho_x	0.093	0.059	0.079	Latent field
23	ui_anc_rho_x	0.037	0.045	0.024	Latent field
24	ui_anc_rho_x	0.056	0.110	0.107	Latent field
25	ui_anc_rho_x	0.031	0.070	0.088	Latent field
26	ui_anc_rho_x	0.032	0.046	0.091	Latent field
27	ui_anc_rho_x	0.065	0.065	0.046	Latent field
28	ui_anc_rho_x	0.026	0.036	0.019	Latent field
29	ui_anc_rho_x	0.087	0.050	0.044	Latent field
30	ui_anc_rho_x	0.018	0.033	0.029	Latent field
31	ui_anc_rho_x	0.028	0.030	0.062	Latent field
32	ui_anc_rho_x	0.031	0.058	0.048	Latent field
1	ui_anc_alpha_x	0.048	0.076	0.065	Latent field
2	ui_anc_alpha_x	0.035	0.076	0.073	Latent field
3	ui_anc_alpha_x	0.061	0.089	0.085	Latent field
4	ui_anc_alpha_x	0.071	0.070	0.061	Latent field
5	ui_anc_alpha_x	0.029	0.049	0.043	Latent field
6	ui_anc_alpha_x	0.042	0.020	0.051	Latent field
7	ui_anc_alpha_x	0.067	0.051	0.082	Latent field
8	ui_anc_alpha_x	0.057	0.048	0.096	Latent field
9	ui_anc_alpha_x	0.053	0.087	0.099	Latent field
10	ui_anc_alpha_x	0.030	0.047	0.060	Latent field
11	ui_anc_alpha_x	0.041	0.026	0.048	Latent field
12	ui_anc_alpha_x	0.043	0.063	0.075	Latent field
13	ui_anc_alpha_x	0.054	0.062	0.074	Latent field
14	ui_anc_alpha_x	0.059	0.073	0.072	Latent field
15	ui_anc_alpha_x	0.048	0.077	0.067	Latent field
16	ui_anc_alpha_x	0.103	0.091	0.046	Latent field
17	ui_anc_alpha_x	0.055	0.037	0.046	Latent field
18	ui_anc_alpha_x	0.133	0.145	0.041	Latent field
19	ui_anc_alpha_x	0.059	0.080	0.069	Latent field
20	ui_anc_alpha_x	0.093	0.058	0.098	Latent field
21	ui_anc_alpha_x	0.078	0.058	0.064	Latent field
22	ui_anc_alpha_x	0.079	0.060	0.036	Latent field
23	ui_anc_alpha_x	0.044	0.042	0.072	Latent field
24	ui_anc_alpha_x	0.062	0.052	0.061	Latent field
25	ui_anc_alpha_x	0.116	0.079	0.073	Latent field
26	ui_anc_alpha_x	0.076	0.045	0.061	Latent field
27	ui_anc_alpha_x	0.071	0.050	0.097	Latent field
28	ui_anc_alpha_x	0.026	0.048	0.039	Latent field
29	ui_anc_alpha_x	0.126	0.122	0.046	Latent field
30	ui_anc_alpha_x	0.070	0.041	0.079	Latent field
31	ui_anc_alpha_x	0.044	0.046	0.059	Latent field
32	ui_anc_alpha_x	0.045	0.052	0.040	Latent field
1	log_sigma_or_gamma	0.524	0.524	0.036	Hyper
1	log_or_gamma	0.035	0.121	0.103	Latent field
2	log_or_gamma	0.033	0.086	0.086	Latent field
3	log_or_gamma	0.049	0.030	0.042	Latent field
4	log_or_gamma	0.053	0.097	0.072	Latent field
5	log_or_gamma	0.060	0.102	0.098	Latent field
6	log_or_gamma	0.060	0.084	0.087	Latent field
7	log_or_gamma	0.036	0.024	0.027	Latent field
8	log_or_gamma	0.058	0.051	0.057	Latent field
9	log_or_gamma	0.018	0.087	0.070	Latent field
10	log_or_gamma	0.064	0.053	0.044	Latent field
11	log_or_gamma	0.035	0.038	0.050	Latent field
12	log_or_gamma	0.029	0.136	0.112	Latent field
13	log_or_gamma	0.041	0.075	0.059	Latent field
14	log_or_gamma	0.035	0.022	0.030	Latent field
15	log_or_gamma	0.026	0.049	0.051	Latent field
16	log_or_gamma	0.052	0.025	0.032	Latent field
17	log_or_gamma	0.036	0.050	0.052	Latent field
18	log_or_gamma	0.034	0.149	0.137	Latent field
19	log_or_gamma	0.045	0.033	0.056	Latent field
20	log_or_gamma	0.030	0.025	0.042	Latent field
21	log_or_gamma	0.072	0.045	0.070	Latent field
22	log_or_gamma	0.024	0.037	0.044	Latent field
23	log_or_gamma	0.043	0.027	0.031	Latent field
24	log_or_gamma	0.066	0.047	0.054	Latent field
25	log_or_gamma	0.048	0.043	0.044	Latent field
26	log_or_gamma	0.052	0.048	0.059	Latent field
27	log_or_gamma	0.025	0.035	0.021	Latent field
28	log_or_gamma	0.077	0.060	0.054	Latent field
29	log_or_gamma	0.042	0.043	0.052	Latent field
30	log_or_gamma	0.036	0.077	0.086	Latent field
31	log_or_gamma	0.053	0.072	0.065	Latent field
32	log_or_gamma	0.043	0.032	0.052	Latent field

This plot could be improved by making the red to green transition smooth, then adding legend for KS difference. The density plots could be filled according to KS difference.

summary_ks_plot <- function(ks_summary, method1, method2) {
  ks_method1 <- paste0("KS(", method1, ", tmbstan)")
  ks_method2 <- paste0("KS(", method2, ", tmbstan)")
  
  densityplot <- ggplot(ks_summary, aes(x = .data[[ks_method1]], y = .data[[ks_method2]])) +
    stat_density_2d(aes(linetype = Type), col = "black") +
    annotate(geom = "polygon", x = c(-Inf, Inf, Inf), y = c(-Inf, Inf, -Inf), fill = multi.utils::cbpalette()[2], alpha = 0.2) +
    annotate(geom = "polygon", x = c(-Inf, Inf, -Inf), y = c(-Inf, Inf, Inf), fill = multi.utils::cbpalette()[6], alpha = 0.2) +
    xlim(0, 0.75) +
    ylim(0, 0.75) +
    geom_abline(intercept = 0, slope = 1, linetype = "dashed", alpha = 0.5) +
    scale_linetype_manual(values = c("dashed", "solid")) +
    # scale_alpha_discrete(range = c(0.5, 0.01)) +
    # guides(alpha = FALSE) +
    guides(linetype = FALSE) +
    labs(x = ks_method1, y = ks_method2) +
    theme_minimal() 

  ks_summary[["KS difference"]] <- ks_summary[[ks_method1]] - ks_summary[[ks_method2]]
  
  ridgeplot <- ggplot(ks_summary, aes(y = Type, x = `KS difference`)) +
    ggridges::geom_density_ridges(alpha = 0.7, fill = NA, aes(linetype = Type)) +
    coord_flip() +
    scale_linetype_manual(values = c("dashed", "solid")) +
    labs(y = "", x = paste0(ks_method1, " - ", ks_method2)) +
    guides(linetype = FALSE) +
    theme_minimal()

  densityplot + ridgeplot
}

summary_ks_plot(ks_summary, "TMB", "aghq")

## Picking joint bandwidth of 0.00943

summary_ks_plot(ks_summary, "TMB", "adam")

## Picking joint bandwidth of 0.00997

summary_ks_plot(ks_summary, "aghq", "adam")

## Picking joint bandwidth of 0.241

(ks_summary_out <- ks_summary %>%
  group_by(Type) %>%
  summarise(
    TMB = mean(`KS(TMB, tmbstan)`),
    aghq = mean(`KS(aghq, tmbstan)`),
    adam = mean(`KS(adam, tmbstan)`)
  ))

Save some things for output:

pdf("ks_summary_tmb_adam.pdf", h = 4, w = 6.25)
summary_ks_plot(ks_summary, "TMB", "adam")

## Picking joint bandwidth of 0.00997

dev.off()

## quartz_off_screen 
##                 2

saveRDS(ks_summary_out, "ks_summary_out.rds")

3.3 Investigation into large KS values

The following parameters have KS values greater than 0.5 in both dimensions.

(big_ks <- ks_summary %>%
  filter(`KS(aghq, tmbstan)` + `KS(TMB, tmbstan)` > 0.5) %>%
  pull(Parameter))

##  [1] "logit_phi_rho_x"      "log_sigma_rho_x"      "logit_phi_rho_xs"     "log_sigma_rho_xs"     "logit_phi_rho_a"      "log_sigma_rho_a"      "logit_phi_rho_as"     "log_sigma_rho_as"    
##  [9] "log_sigma_rho_xa"     "logit_phi_alpha_x"    "log_sigma_alpha_x"    "logit_phi_alpha_xs"   "log_sigma_alpha_xs"   "logit_phi_alpha_a"    "log_sigma_alpha_a"    "logit_phi_alpha_as"  
## [17] "log_sigma_alpha_as"   "log_sigma_alpha_xa"   "OmegaT_raw"           "log_betaT"            "log_sigma_lambda_x"   "log_sigma_ancrho_x"   "log_sigma_ancalpha_x" "log_sigma_or_gamma"

4 MMD

#' To write!

5 PSIS

Suppose we have two sets of samples from the posterior. For each sample we are going to want to evaluate the log-likelihood, so that we can calculate the log-likelihood ratios. We can extract the TMB objective function for the log-likelihood as follows:

tmb$fit$obj$fn()

## [1] NaN

#' To write!